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Abstract 

The  multidimensional  inverse  scattering  problem  for  an 
acoustic  medium  is  considered  within  the  homogeneous 
background  Bom  approximation.  A  constant  density 
acoustic  medium  is  probed  by  a  wideband  plane  wave 
source,  and  the  scattered  held  is  observed  along  a  re¬ 
ceiver  array  located  outside  the  medium.  The  inversion 
problem  is  formulated  as  a  generalized  tomographic 
problem.  It  is  shown  that  the  observed  scattered  field 
can  be  appropriately  filtered  so  as  to  obtain  generalized 
projections  of  the  scattering  potential.  For  a  2-D  exper¬ 
imental  geometry,  these  projections  are  weighted  inte¬ 
grals  of  the  scattering  potential  over  regions  of  parabolic 
support.  The  inversion  problem  is  therefor  similar  to 
that  of  x-ray  tomography,  except  that  instead  of  be¬ 
ing  given  projections  of  the  object  to  be  reconstructed 
along  straight  lines,  projections  along  parabolas  are 
given.  The  inversion  procedure  that  we  propose  is  sim¬ 
ilar  to  the  x-ray  solution,  in  the  sense  that  it  consists 
of  a  backprojection  operation  followed  by  2-D  space  in¬ 
variant  filtering.  A  "Projection-Slice  Theorem"  is  de¬ 
rived  relating  the  generalized  projections  and  the  scat¬ 
tering  potential  in  the  Fourier  transform  domain. 

1.  Introduction 

In  inverse  scattering  problems,  the  objective  is  to 
reconstruct  certain  physical  properties  of  a  medium 
from  scattering  experiments.  In  general,  there  is  an 
array  of  sources  and  an  array  of  receivers  located  out¬ 
side  the  medium.  There  are  three  general  approaches 
to  the  inversion  problem:  Iterative  inversion,  exact  di¬ 
rect  inversion  and  approximate  direct  inversion.  The 
iterative  inversion  (also  called  generalized  inversion) 
approach  attempts  to  solve  a  very  large  scale  least- 
squares  problem,  where  at  every  stage  an  estimate  of 
the  medium  parameters  is  used  to  solve  the  forward 
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scattering  problem  for  the  corresponding  wave  field  at 
the  receiver  locations.  Depending  on  the  difference  be¬ 
tween  the  observed  and  the  computed  scattered  waves, 
a  new  estimate  of  the  medium  model  is  obtained  and 
the  next  iteration  is  performed.  These  methods  axe 
very  time  consuming,  and  the  convergence  of  the  cur¬ 
rent  methods  depends  on  the  accuracy  of  the  initial 
estimate. 

The  objective  of  exact  direct  inversion  methods  is 
to  reconstruct  the  medium  properties  exactly,  with  no 
iterations  involved.  These  methods  require  large  num¬ 
bers  of  sources  and  receivers  with  particular  observa¬ 
tion  geometries  which  limit  their  applicability  to  prac¬ 
tical  problems.  From  a  theoretical  point  of  view,  one¬ 
dimensional  exact  direct  inverse  scattering  methods 
have  reached  an  advanced  level  of  development  (see  [lj 
for  an  overview),  whereas  their  extension  to  higher  di¬ 
mensions  has  proved  to  be  difficult. 

Consequently,  the  logical  approach  to  the  practi¬ 
cal  solution  of  the  multidimensional  inverse  scatter¬ 
ing  problems  is  the  use  of  approximate  direct  inversion 
methods.  The  differential  equation  for  wave  propaga¬ 
tion  in  a  medium  can  be  transformed  into  the  so-called 
Lippmann-Schwinger  integral  equation  [12j.  For  ex¬ 
ample,  for  an  acoustic  medium  with  constant  density, 
this  equation  expresses  in  integral  form  the  scattered 
field  inside  the  medium  in  terms  of  the  propagation 
velocity  profile  and  the  total  pressure  field  inside  the 
medium.  The  Bom  approximation  consists  of  approxi¬ 
mating  the  total  field  inside  the  integral  representation 
by  the  incident  field.  Therefore,  this  approximation  as¬ 
sumes  that  the  scattered  field  is  small  compared  to  the 
incident  field;  or  equivalently  that  the  perturbations  of 
the  velocities  are  small  with  respect  to  the  background 
velocity  profile  which  is  used  to  compute  the  incident 
field.  Another  way  of  interpreting  the  Bom  approxi¬ 
mation  is  to  view  it  as  a  linearization  technique,  where 
the  relation  between  the  scattered  field  and  the  object 
profile  that  we  want  to  reconstruct  is  linear.  Physi¬ 
cally,  the  Bom  approximation  takes  into  account  only 
the  singly  scattered  waves;  multiply  scattered  waves 
are  considered  as  noise.  Note  that,  depending  on  the 
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method  used  to  compute  the  background  field,  the  mul¬ 
tiples  due  to  the  background  model  may  be  included  in 
the  scattered  field.  The  multiples  due  to  the  residual 
velocities  are  neglected. 

In  this  paper,  the  inverse  scattering  problem  for  an 
acoustic  medium  is  considered  within  the  homogeneous 
background  Bom  approximation.  A  constant  density 
acoustic  medium  is  probed  by  a  wide-band  plane  wave 
source,  and  the  scattered  field  is  observed  along  a  re¬ 
ceiver  array  located  outside  the  medium.  The  objec¬ 
tive  is  to  reconstruct  the  scattering  potential,  which  Is 
a  function  of  the  propagation  velocity  inhomogeneities 
in  the  medium.  The  monochromatic  plane-wave  source 
inverse  scattering  problem  has  been  studied  under  the 
name  of  diffraction  tomography  by  several  people,  in¬ 
cluding  Mueller,  Kaveh,  Devaney  and  Beylkin  [3,4,7], 
Roberts  and  Kak  investigated  the  reflection  mode  prob¬ 
lem  with  broadband  illumination  [ll],  whereas  Esmer- 
soy  and  Levy  presented  a  solution  in  terms  of  an  ex¬ 
trapolated  field  [5j. 

In  the  present  paper,  the  key  observation  is  that, 
the  time  traces  observed  at  the  receivers  can  be  ap¬ 
propriately  filtered  so  as  to  obtain  generalized  projec¬ 
tions  of  the  scattering  potential.  For  a  two-dimensional 
experimental  geometry,  these  projections  are  weighted 
integrals  of  the  scattering  potential  over  regions  of  pa¬ 
rabolic  support,  whereas  they  become  surface  integrals 
over  paraboloids  for  the  three-dimensional  case.  Thus 
the  inverse  scattering  problem  is  now  posed  as  a  gen¬ 
eralized  tomographic  or  integral  geometric  problem. 

The  straight-line  tomography  problem,  which  arises 
in  x-ray  tomography,  was  first  solved  by  Radon  [10] ;  see 
Deans  [2]  for  a  full  treatment.  Fawcett  [6]  formulated 
the  zero-offset  Bom  inversion  problem  as  a  generalized 
tomographic  problem,  where  the  objective  is  to  recon¬ 
struct  a  function  from  its  projections  along  circles  or 
spheres. 

The  reconstruction  problem  for  parabolical  projec¬ 
tions  can  be  formulated  in  a  way  similar  to  the  problem 
of  x-ray  tomography.  The  solution  can  be  expressed  as 
a  backprojection  operation  where  we  sum  the  contribu¬ 
tions  of  all  projections  passing  through  a  given  point  in 
space,  followed  by  a  two  or  three  dimensional  filtering 
operation. 

A  different,  physically-oriented  interpretation  of  the 
backprojection  operation  appearing  in  our  reconstruc¬ 
tion  technique  can  be  developed  by  showing  that  it 
corresponds  to  first  backpropagating  the  observed  fil¬ 
tered  time  traces,  and  then  correlating  the  backprop- 
agated  field  with  the  incident  probing  wave  field.  The 
details  of  this  interpretation  will  be  reported  in  [8].  A 
"Projection-Slice  Theorem"  is  also  derived  relating  the 
generalized  projections  and  the  scattering  potential  in 
the  Fourier  transform  domain. 


Fig.  I  The  experimental  geometry. 


2.  Problem  Description 

In  this  paper  we  will  treat  the  two-dimensional  case; 
a  more  detailed  account  including  the  three-dimensional 
case  will  be  reported  elsewhere  [8j.  Consider  the  scat¬ 
tering  experiment  described  in  Fig.  1.  A  constant  den¬ 
sity  2-D  acoustic  medium  is  probed  by  a  wide-band 
plane  wave  and  the  scattered  field  is  observed  along 
a  straight-line  receiver  array.  The  Fourier  transform 
P(z,ui)  of  the  pressure  field  at  position  x=  (x,y)  sat¬ 
isfies 

[V*  4-  k2n2(x}] P(z,  ui)  =  0,  (l) 

where  k  =  w/c0  is  the  wavenumber,  rt(x)  =  c0/c(x)  is 
the  refractive  index  of  the  medium,  c(x)  is  the  prop¬ 
agation  velocity  at  point  x  and  Co  is  the  propagation 
velocity  of  the  background  medium.  We  assume  that 
n(x)  does  not  deviate  significantly  from  the  background 
index  1,  so  that 

n2(z)  =  1  +  U(x),  (2) 

where  the  scattering  potential  I/(x)  is  small.  We  also 
assume  that  U(x)  has  a  bounded  support  V,  which 
is  completely  located  on  the  same  side  of  the  receiver 
array,  as  shown  in  Fig.  1.  Decomposing  the  total  field 
P  into  the  incident  field  Pq  and  the  scattered  field 

P(x,  w)  =  P0(x,  w)  +  P,(x,  w),  (3) 

and  noting  that  the  incident  field  Po  satisfies  the  Helm¬ 
holtz  equation  associated  with  the  background  medium 
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(V2  +  *Vo(®.«)-0,  (4) 

we  find  that  the  scattered  field  P,  satisfies  the  equation 
(V2  +  k2)P,(xj  w)  =  -fc2tr(x)P(*,u;).  (5) 

The  solution  of  (5)  is  given  by  the  Lippmann-Schwinger 
equation  [12j 

P,(Xju)=k2J  dxfZ7(P)G0(£xf,w)P(P,w),  (6) 

where  Gq(x,P,u;)  is  the  Green’s  function  associated 
with  a  point  source  in  a  homogeneous  medium: 

(V2 +  fc2)G0(£,  2?,  «)=-*($- sO-  (7) 

Equation  (6)  demonstrates  the  nonlinear  relation 
that  exists  between  the  potential  U (x)  and  the  pressure 
field  P(z,ui).  To  linearize  this  equation  we  adopt  the 
Bom  approximation,  whereby  we  assume  P,(x,  ui)  <?C 
Po(xjUi).  Hence  the  Lippmann-Schwinger  equation  be¬ 
comes 

•f*(£>w)  «  k2  J dsfCT(P)G0(£jS!j«)-f,o(sf7w)-  (8) 

For  the  problem  under  consideration,  the  Green’s  func¬ 
tion  and  the  incident  wave  are  given  as 

G0(x,s',ci;)  =  iH«(%-S'|),  (9) 

P0(x»  =  e“K  (10) 

where  £  is  the  unit  vector  which  indicates  the  single  of 
incidence  of  the  plane- wave  source,  and  Pg^(-)  indi¬ 
cates  the  Hankei  function  of  order  zero  and  type  one. 
Therefore,  the  scattered  field  at  a  receiver  point 

£  =  (pcosfli  +  £sin<£,  p sin 4>  —  £cos<£) 

is  given  by 

*.(£.«)  =  V/ (11) 

within  the  Bom  approximation.  Hence,  we  have 

3  >(€,*)-  (12) 

Define  the  inverse  Fourier  transform  of  P(£,k)  with 
respect  to  k  as  g(  £,  r) : 

g(z,r)  =  f^dkF{Z,k)e-ik'.  (13) 

Taking  account  of  the  fact  that  (see  [9j,  p.  731) 


Fig.  2  The  generalized  parabolic  projections. 


where  l(-)  i3  the  unit  step  function,  therefore  we  find 
that 


9{Z,r)  =  /  #l7(af) 


l(r-£-P-[sf- £|) 
v/(r_£.2/)z_|x'_^|2’ 


(15) 


Thi3  identity  expresses  y(£,  r)  as  a  weighted  integral 
of  the  scattering  potential  U[x)  where  the  weighting 
function  is  nonzero  in  a  region  with  parabolic  support. 
The  parabola  satisfies  the  equation  r  =  £  •  £  —  \s.  —  £| 
where  r,  £  and  £  are  given  and  x,  varies,  as  shown  in 
Fig.  2.  The  weighting  function  becomes  infinite  for 
values  of  x  along  this  parabola,  so  that  the  largest  con¬ 
tribution  to  the  integral  is  made  by  the  values  of  U (r) 
which  lie  along  the  parabola.  In  some  sense,  ?(£,r)  is 
a  projection  of  U{xj  with  respect  to  a  function  whose 
singularities  are  algebraic  and  located  along  a  parabola. 

It  is  then  interesting  to  note  that  the  projections 
ff(£j r)  can  be  obtained  directly  from  the  observed  scat¬ 
tered  field  P,(£,  f)  in  time  domain  also:  From  (12)  and 

(13), 


Pt/cq  [T 

g(Z,r)  =  -2*c0  dr  dsP,(£,s).  (16) 

Jo  Jo 


Thus,  the  projections  p(£,r)  are  proportional  to  the 
scattered  field  twice  integrated. 

Therefore,  given  the  observed  scattered  field  P,(£,  t) 
for  £  along  the  line  array  depicted  in  Fig.  2,  and  for 
— oo  <  t  <  oo,  it  is  a  straightforward  procedure  to  find 
the  generalized  projection  p(£,r)  for  — oo  <  £  <  oo 
and  0  <  r  <  oo,  and  in  the  following  it  will  be  assumed 
that  these  projections  constitute  the  data  that  is  given 
by  the  scattering  experiment.  From  this  point  of  view, 
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the  inverse  scattering  problem  can  be  formulated  as 
follows:  Given  the  generalized  projections 

{ff($,  i* )  '  — oo  <$<oo,  0  <  r  <  oo}, 

we  want  to  reconstruct  the  scattering  potential  U(xj. 

In  some  sense,  the  above  reconstruction  problem  is 
similar  to  the  problem  of  x-ray  tomography,  where  we 
are  given  the  projections 

g{r,9)  =  J  dJU{£)6(r-h£) 

of  Z7(z)  along  straight  lines.  It  is  also  analog  to  the 
inverse  scattering  problem  considered  by  Fawcett  [6], 
where  it  was  shown  that  the  so-called  zero-offset  inverse 
scattering  geometry  could  be  reduced  to  the  solution  of 
a  generalized  tomographic  problem  where  the  objective 
is  to  reconstruct  a  function  U(x)  from  its  projections 
along  circles  of  arbitrary  radii  centered  along  a  straight 
line.  There  exists  however  an  important  difference  be¬ 
tween  the  problem  of  x-ray  tomography,  or  the  prob¬ 
lem  considered  by  Fawcett  [6|,  and  the  problem  that 
we  examine  here.  In  x-ray  tomography,  the  projections 
of  U(xj  are  taken  with  respect  to  a  weighting  function 
which  is  the  distribution  <5(r— £-x),  whereas  the  weight¬ 
ing  function  appearing  in  (15)  is  algebraic.  Thus,  in  2- 
D  the  generalized  parabolic  projection  gr($,r)  is  not  an 
integral  along  a  curve  but  an  integral  that  has  nonzero 
weighting  over  a  whole  region  of  the  plane  (the  inside 
of  the  parabola  depicted  in  Fig.  2).  This  is  due  to 
the  fact  that  the  2-D  Green’s  function  does  not  have 
an  impulsive  waveform,  but  has  a  tail,  as  indicated  by 
equation  (9). 


3.  The  Backprojection  Operation 


Like  the  x-ray  tomography  inversion  procedure,  or 
the  method  proposed  by  Fawcett  [6]  for  the  case  of 
circular  projections,  the  first  step  of  our 'inversion  pro¬ 
cedure  is  to  perform  a  backprojection  operation  on  the 
projections  p($,  r).  At  a  given  point  x,  this  operation 
sums  the  contributions  of  the  projections  ff($,  r)  which 
correspond  to  parabolic  regions  that  include  the  point 
x.  In  this  summation,  the  projection  <;($,  r)  is  weighted 
in  proportion  to  the  amount  U{xj  has  contributed  to 
it  according  to  the  forward  scattering  equation  (15). 
By  performing  this  backprojection  operation  for  every 
point  in  the  plane,  this  gives  a  backprojection  approx¬ 
imation,  Ub(x)  to  U(x).  It  is  given  by 


*s(x) .  r  d$ 

J-°°  Jo  y/{r  -9 ■£)*  -|£- $|2 

(17) 

Our  first  objective  is  to  find  a  frequency  domain  re¬ 


lationship  between  <?($,  r)  and  UB(x}.  It  can  be  shown 
that  [8] 


/ ;r-*«  , 
•>  yj  (r  -  £•  z)2  -  |x- $|2 


ITT 


.  ( k2  -  k2)  cos  (9  +  <f>)  4-  2 kxkv  sin  ( 9  +  <f>) 
eX^>  ^  tp  2  (k*  cos  9  +  kv  sin  9) 

{k\  -  kj)  sin  [9  +  <j>)  —  2k,kv  cos  {9  +  <j>) 
■  exp  <  -»$ - - - 


£|2  kxcos9  +  kvsm9 

} 


2(k*  cos  9  +  kv  sin  9) 
k 2 

2(kx  cos  9  H-  ktf  sin  9) 
where  fc  =  (k*,k„)  and  k  =  |fc|.  Therefore 


•  exp  \  —  tr 


}■ 


(18) 


f  .  (k*  —  k^  cos  {9  +  4))+  2 kxkv  sin  {9  +  <f>)  | 
^  (  ***  2(k*  cos  9  +  kv  sin  9)  J 


2(k^  cos  9  +  k9  sin  9) 

-  f  _  {kl  -  kl)  sin  {9  +  <t>)  —  2kxkv  cos  {9  -+•  <f>) 
?  2(k„cos0  +  k9sin0) 

k2 


K  = 


i9))' 


(19) 


2(ks  cos  9  +  fcjsinl 
where  j(k£,kr)  is  the  2-D  Fourier  transform  of  p($,  r). 


4.  The  Relationship  Between  Ug{k) 
and  U{k) 

We  first  take  the  Fourier  transform  of  p($,r)  with 
respect  to  r: 

ff($.*r)  =  /  dre~ikr'g{£,r)  (20) 

Jo 

=  F'{M 

=  -y/ dz'C/(z')eii^'4J,(k|x'  -  $|), 

from  (12).  Now  taking  the  Fourier  transform  with  re¬ 
spect  to  $,  we  obtain  [9] 

|  (k?,kr)  =  -  t7^3f1^r) 
yjk;  -  ^ 

•U{kx  =  kr  cos  9  +  k^  sin  <j>  +  Tfsgn (kr) yjkj  —  ki  cos  c/>, 
ky  —  kr  sin  9  —  fcf  cos  <j>  +  7sgn(kr)  \Jkf— k|  sin  4>) , 
for  |fc€|  <  |k,|, 

§(fc£,kr)  =  0,  for  |k£|  >  \kr\,  (21) 

and  where 
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XJ(]A  —  ~  K)  coa  (9  +  &)  +  2kTkv  sin  (8  4-  <£)] 

~  2?r(A:I  cos  8  +  ky  sin  8) 

—  k2)  cos  (8  -+•  <j>)  +  2kxky  sin  ( 8  +  <f>)\ 


•exp 


f  •  (*5 
l"‘p- 


2  (A*  cos  5  +  ky  sin  8) 


*  .  £  _  (fcg  —  fey)  sin  (fl  +  <ft)  —  2kxky  cos  (8  +  <j>) 

^  4  2(ArIcos0  +  fc„sin0)  ’ 


kr  = 


k2 


2(kz  cos  8  -rky  sin  i 
Combining  (19)  and  (23)  gives 

frua- 


,*€<7.  ■  (23) 


2^3  [(**  ~  cos  (a  +  p)  +  2*x*w  sin  (5  +  <£)]  Dfl(i), 

.  (24) 

for  k£C.  Therefore,  U (&)  can  be  obtained  from  Ub{M} 
by  a  single  2-D  filtering  operation.  This  is  similar  to 
the  "filter  of  backp  rejections”  method  used  in  straight- 
line  tomography  (see  [2],  Section  6.5). 

By  using  the  generalized  parabolical  projections 
g{£,  r)  obtained  from  a  single  straight-line  array,  we  ob¬ 
tain  the  coverage  of  U(k)  over  a  90°  cone.  To  obtain  a 
more  complete  coverage,  we  can  use  additional  receiver 
arrays  or  perform  several  experiments  with  plane-waves 
incident  from  several  directions. 


+1  if  ( x cos<j>  +  y  sin<£  —  p)  >  0  for  all  x, y  €  V, 
—1  if  (xcos<£  +  ysin<£  —  p)  <  0  for  all  x,y  €  V. 

For  normal  incidence  (8  =  <j>)  and  fixed  fc,,  this  result 
corresponds  to  the  formula  obtained  by  Devaney  in  the 
context  of  diffraction  tomography  [3].  This  relationship 
is  similar  to  the  Projection  Slice  Theorem  of  straight- 
line  tomography  (see  [2j,  Section  6.2).  It  relates  the  1-D 
Fourier  transform  of  £(£,  kr)  with  respect  to  £,  and  the 
2-D  Fourier  transform  of  U(x}.  For  a  fixed  kr,  g(k kr) 
gives  U(k)  along  two  semicircles  of  radius  kr  centered 
at  (krcosd,kram8)  and  {—k,  cos  8,  —kr  sin  8)  (see  Fig. 
3).  By  letting  kr  vary,  we  see  that  U(k)  is  only  known 
in  the  cone  C ,  which  is  defined  as 

C  = 

{*„  ky  :  i  (e  +  <j>  -  |)  <  tan-1  jjj  <  \  (*  +  <j>  +  |) 

(22) 

for  7  —  +1.  For  7  =  —1,  C  is  the  complement  of  the 
above. 

The  inverse  formula  of  (21)  is 


5.  Numerical  Example 

We  present  computer  simulation  results  for  simula¬ 
tion.  Fig.  4  shows  the  scattering  potential  to  be  recon¬ 
structed.  It  represents  a  ±5%  variation  in  the  propaga¬ 
tion  velocity  with  respect  to  the  constant  background 
velocity.  In  the  experiments,  the  probing  plane  wave 
is  incident  perpendicularly  to  Side  A  and  the  receivers 
were  placed  on  all  four  sides.  The  scattered  waves  were 
computed  by  a  finite-difference  algorithm.  The  source 
wavelet  used  was  a  Black man-Harris  window  and  the 
diameter  of  the  object  was  twice  the  dominant  wave¬ 
length.  Fig.  5  shows  the  reconstruction  obtained  using 
only  the  receivers  on  Side  A,  while  Figs.  6  and  7  show 
the  same  for  the  cases  where  the  receivers  are  only  lo¬ 
cated  on  Sides  B  and  C  respectively.  Fig.  8  depicts  the 
reconstruction  using  all  the  receivers. 

6.  Conclusion 

In  this  paper  we  have  considered  the  direct  velocity  in¬ 
version  problem  for  a  constant-density  acoustic  medium 
probed  by  a  single  wide-band  plane  wave.  The  prob¬ 
lem  was  posed  as  a  generalized  tomographic  problem, 
where  weighted  integrals  of  the  scattering  potential 
,U(s}  over  parabolic  regions  Me  considered  as  data. 
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Drawing  analogy  to  x-ray  tomography,  a  backprojec- 
tion  operator,  Ub(x)  was  defined,  which  can  be  viewed 
as  a  "migration”  approximation  to  U(x).  Ub(x)  was 
related  to  the  parabolical  projections  g(£,  r)  in  the  2-D 
Fourier  transform  domain.  The  parabolical  projections 
were  also  related  to  U (2)  in  the  frequency  domain,  thus 
deriving  a  "Projection  Slice  Theorem". 


Fig.  6  Reconstruction  using  the  receivers  on  Side  B. 
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